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'/Describes  the  mathematics  and  the  computer 
program  for  solving  the  problem  of  optimum 
bonus  management  formulated  ir^Patricia 
Munch* 'A  Dynamic  Model  for  Optimum  Bonus 
Management',  R-T940*frRPA7  danuary~1977i  The 
problem  of  optimum  bonus  management  is  treated 
as  a discrete  linear  control  system  with  a 
quadratic  cost  function  and  solved  by  using 
Pontryagin's  discrete  maximum  principle. 

The  state  of  the  system  at  discrete  time 
is  a vector  of  numbers  of  men  in  each  of  a 
set  of  year  groups.  The  system  evolves 
linearly  in  time  under  linear  controls  that 
are  the  bonuses  paid  to  the  men  in  a pre- 
scribed subset  of  the  year  groups.  The 
program  solves  for  the  sequence  of  bonus 
values  that  drive  a given  inital  state  of 
year  groups  to  a prescribed  final  state  and 
minimize  a sum  of  quadratic  bonus  and  penalty 
costs.  The  program,  which  consists  of  a main 
program  and  22  subroutines,  is  written  in 
FORTRAN  IV  double  precision  for  the  IBM 
370-158.  (Author) 
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PREFACE 


This  report  was  prepared  under  Rand's  DoD  Training  and  Manpower 
Program,  sponsored  by  the  Cybernetics  Technology  Office  of  the  Defense 
Advanced  Research  Projects  Agency  (ARPA).  The  purpose  of  this  re- 
search program  is  to  develop  broad  strategies  and  specific  solutions 
for  dealing  with  present  and  future  military  manpower  problems,  includ- 
ing the  development  of  new  research  methodologies  for  exar.ining  broad 
classes  of  manpower  problems,  as  well  as  specific  problem-oriented  re- 
search. It  is  hoped  that  this  research  program  will  provide  analysis 
of  current  and  future  manpower  issues  and  contribute  to  a better  under- 
standing of  the  manpower  problems  confronting  the  Department  of  Defense. 

This  report  is  a companion  to  Patricia  Munch,  A Dynamic  Model  for / 
i'l'tirKtn  BcnuS  tiina-jemcnt , R-1940-DoD/ARPA,  January  1977,  which  deals  ' 
with  the  institutional  background  of  reenlistment  bonus  management  and 
its  formulation  both  as  an  economic  problem  and  as  a mathematical  prob- 
lem of  optimal  control. 

The  present  report  describes  the  mathematical  foundations,  oper- 
ating instructions,  and  structure  and  function  of  the  Dynamic  Bonus 
Management  Model  Program  (DBMMP)  and  its  subroutines.  Its  purpose  is 
to  give  a clear,  detailed  derivation  of  the  formulas  used  by  the  pro- 
gram to  determine  the  optimal  bonus  policy.  It  provides  complete, 
clear  DBMMP  input  parameter  definitions  and  input  parameter  data  card 
format  specifications  that  will  permit  the  user  to  run  the  program 
with  ease.  The  report  also  supplies  the  programmer  with  sufficient 
knowledge  of  program  function,  architecture,  and  logical  flow  to  modify 
the  program,  as  required,  for  greater  generality. 
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SUMMARY 


This  report  provides  a detailed  mathematical  derivation  of  the 
algorithm  for  solving  the  optimal  bonus  management  problem  that  was 
posed  in  Patricia  Munch,  A M Id  j'  ■>•  >ftitnur  Mana<ffm>  i:t, 

R- 1 940-I)ol)/ARPA , January  1977.  In  addition  to  the  mathematical  de- 
velopment, the  report  also  describes  the  computer  program  that  incor- 
porates the  solution  method  and  defines  the  economic  setting  of  the 
problem . 

The  problem  of  optimum  bonus  management  is  treated  as  a discrete 
linear  control  system  in  which  the  optimal  control  is  determined  from 
Pontryagin's  discrete  maximum  principle. 

In  this  formulation  the  state  of  the  system  at  time  t is  a vector, 
each  of  whose  components  is  the  number  of  men  in  each  of  a number  of 
year  groups.  The  system  evolves  linearly  in  discrete  time  steps  under 
the  linear  action  of  a set  of  control  variables.  The  control  variables 
are  bonuses  that  are  paid  to  the  men  in  certain  year  groups. 

The  optimization  problem  involves  the  determination  of  the  se- 
quence  of  bonuses  that  will  drive  prescribed  initial  state  of  year 
group  distribution  to  a desired  terminal  state  in  a given  number  of 
years  while  minimizing  a cost  function.  The  cost  function  is  a cumu- 
lative sum  of  two  quadratic  forms,  penalty  cost  and  bonus  cost.  The 
penalty  cost  quadratic  form  represents  the  costs  accruing  to  shortages 
or  overages  in  the  year  group  numbers  relative  to  the  corresponding 
desired  terminal  state  values. 


Both  penalty  and  cost  quadratic  term  components  are  functions  of 
a set  of  economic  model  parameters  that  are  defined  in  the  introductory 
section.  The  introduction  also  indicates  which  parameter  values  are 
program  constants  and  which  must  be  supplied  by  the  program  reader. 

The  functional  relation  between  the  economic  model  parameters  and  the 
components  of  the  quadratic  forms  are  developed  in  App.  B. 

The  computer  program  Incorporating  the  solution  method  is  written 
in  FORTRAN  IV  double  precision  for  the  IBM  'J 70-1 ‘>8.  It  contains  a main 


program  and  22  subroutines. 


fRBCEDINQ  PAGI  BLAMUWOT  RIMED 


Although  the  description  in  the  sequel  is  oriented  toward  tile 
solution  of  the  general  bonus  management  optimal  control  problem,  the 
program  is  essciil  ially  a "one-shot"  research  task  and  has  been  "hard 
coded"  for  a system  of  eight  year  groups  with  bonuses  for  men  going 
into  the  first  and  fifth  year  groups.  This  limitation  should  he 
removable  without  a major  reprogramming  effort. 
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I.  INTRODUCTION 


This  report  describes  the  mathematical  solution  method  of  the 
computer  program  that  puts  into  operation  the  model  described  in 
F.  I’.  Munch,  A Dynamic  Model  for  Optimal  Bonus  Management , The  Rand 
Corporation,  R-l  940-I)oI)/ARPA,  January  1977.  The  purpose  of  this  model- 
ing effort  is  to  structure  the  problem  of  bonus  management  as  an  eco- 
nomic problem,  incorporating  into  a single  decisionmaking  framework 
the  several  criteria  currently  used  in  setting  bonuses.  The  dynamic 
adjustment  model  determines  the  optimum  level  of  bonuses  for  several 
year  groups  in  a military  specialty,  at  successive  points  in  the  plan- 
ning horizon,  subject  to  the  constraint  that  the  deviations  of  the 
actual  f rom  the  desired  manpower  inventory  be  eliminated  over  the 
assigned  period.  The  optimum  structure  of  bonuses  minimizes  the  sum 
of  two  sources  of  cost — bonus  cost  incurred  to  reduce  shortages  and 
penalty  cost  assigned  to  shortages  and  overages.  Penalty  cost  is  a 
quadratic  form  in  the  deviations  of  actual  from  desired  manning  levels, 
where  the  quadratic  coefficients  are  derived  from  the  production  and 
demand  functions  for  specialty  output.  These  functions  aggregate  men 
of  differing  skill  levels  into  a homogeneous  unit  of  "output"  (man- 
years)  and  assign  a dollar  value  to  this  output.  Bonus  cost  is 
derived  from  the  supply  function,  which  relates  number  of  men  to  the 
level  of  bonus  paid.  Those  concepts  currently  used  as  rules  of  thumb 
by  bonus  managers,  such  as  criticality  of  a specialty,  substitutability 
between  men  in  different  year  groups,  and  responsiveness  of  supply  to 
a bonus,  are  characterized  by  parameters  in  the  component  functions  of 
the  model.  The  values  of  these  parameters  are  to  be  supplied  by  the 
user.  (See  Table  1.)  The  component  functions  of  the  model  are  re- 
stated below  to  facilitate  specification  of  the  values  of  the  input 
parameters  by  the  user.  Section  II  contains  a statement  of  the  mathe- 
matical problem  and  its  solution. 
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ill  I' IN  I I ION  01  VAK  I Alil.KS 


To  Bf 

Variable  Definition  Supplied  by  User 

= specialty  "output"  in  year  t,  a weighted 
sum  of  number  of  men  in  each  term  of 
serv  ice . '* 


i , l 


X . 
1 


, T 


number  of  man  years  in  ith  term  of  ser- 
v i ce  i n year  t . 


= desired  number  of  man  years  in  ith  term 
of  service. 


A.  = no-bonus  supply  of  lateral  entrants  to 

the  specialty.  In  the  current  model 
only  A j is  non-zero. 


= number  of  men  in  the  jth  year  group 
(YOS)  in  ith  term  of  service  in  year  l. 


,T 


= desired  number  of  men  in  the  jth  year 
group  (YOS)  in  year  T. 


5. 


i 


P 


u 


distribution  parameter, 
tive  productivity  of  men 
terms  of  service.  To 


values,  note  that 

Marginal  product  Xj 
Marginal  product  X 

and  ) 6 . = 1 . 

, i 


reflecting  re  la- 
in different 
select  appropriate 


= substitution  parameter,  defined  as 

follows:  o = 1/(1  + p),  where  a = elas- 
ticity of  substitution  between  X)  and  X> • 
Since  0 5 0 £ + °°,  admissible  values  of 
p run  f rotn  -1  to  + «». 


= returns  to  scale  parameter  p 1 implies 
an  equiproport ionate  increase  in  all  in- 
puts by  x%  leads  to  an  increase  in  out- 
put x/,". 


i 


V.i  r i ab  1 *■ 


To  Pe 

Supp  I i od  by  User 
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T.'ihle  1 --Uont  inueil 


Def in i t ion 


It . 

= no-bonus 

continuation  rate. 

i = 1,2  ... 

) 

7— i .e.  , 

percentage  of  men 

in  YOS  that 

emit  i nue 

into  YOS.  , if  no 

1 + 1 

bonus  is  paid 

B . 

It 

bonus  paid  to  YOS.  in  vear 

1 

t . 

B. 

1 

= slope  of 

supply  function. 

P 

= price  per 

unit  of  specialty 

output . I f 

the  dynamic  adjustment  model  is  run  in 
conjunction  with  a steady-state  model,  P 
is  set  equal  to  the  value  of  the  Lagrange 
multiplier  obtained  by  solving  for  the 
cost-minimizing  input  mix,  subject  to 
producing  the  target  output  level,  7.  . 

If  the  dynamic  model  is  run  without  a 
steady-state  model,  P may  be  assigned  by 
setting  the  target  inventory  and  one  wage 
level,  say  W2 , at  arbitrary  but  reason- 
able levels  and  solving  from  the  first 
order  condition  of  a steady-state  model: 


l/e  = elasticity  of  price  with  respect  to  quan- 

tity. If  a deviation  of  specialty  output 
from  the  target  level  of  xZ  results  in  a 
change  of  more  (less)  than  xl , £ should 
be  set  less  (greater)  than  unity.  £ re- 
flects the  criticality  of  a specialty. 

6 = a scale  parameter.1 


aFor  a discussion  and  economic  interpretation  of  this  production 
function  see  K.  Sato,  "A  Two-Level  Constant-Elast ic ity-of-Subst itution 
Production  Function,"  Reoieti  of  Economic  dtudiosi,  Vol.  34,  April  1967. 

^For  an  alternative  specification,  allowing  differences  in  produc- 
tivity of  men  in  different  year  groups  within  a term  of  service,  see 
F.  P.  Munch,  A Dyndinic  Mod<  ! for  Ik'nuts  Miinajornont. , R- 1 940-l)ol)/ 

AKPA,  January  1977. 

6 is  denoted  <x  in  Munch  (1977).  The  notation  is  changed  here  to 
avoid  confusion  with  the  continuation  rate  parameters. 


II.  THE  OPTIMAL  CONTROL  PROBLEM 


DEFINITION  OF  SYMBOLS 

I. * = | \..)  t»  t Lf  1 ' (’  denoting  transpose) 


is  the  reduced  state  of  the  system  in  year  t. 


»;  ■ 


B2.t ■ 


B , , . . m- 1 . 
1+4 [--£—] , t 


where  |x|  is  the  largest  integer  £ x,  is  the  vector  of  bonus  control 
variables  for  year  n. 


0 

0 

0 

0 


0 

0 

0 


m-l*m-l 


is  the  reduced  state  transition  matrix.  a„ , a_,  ....  a , are  no- 

i i m- 1 

bonus  continuation  rates. 


Although  there  are  m year  groups,  the  dimension  of  the  system 
can  be  reduced  by  one  because  the  first  year  group  number,  Lj  { is  a 
linear  function  of  the  control  variable  B : 

A « t 
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< i B, 


o 

o 


o o . . . o 

0 0 . . . 0 

()  urb,  . . . 

1 ) 


0 


1+' 


0 


m- I * I + | 


111“ 

4 


is  ilu'  control  variables  matrix  multiplier. 


t 


0 J 

m-  1 y 1 


CK, 


is  constant  no-bonus  number  of 
enlistees  into  the  first  YOS) 


is  the  reduced  constant  first  year  group  input. 


0 

0 

• • • 

0 

() 

0 

0 

0 

0 

0 

0 

B5 

0 

0 

. . . 

• > a 

0 

0 

. . . 

• • • 

Bl«| 1 

0 

. . . 

. . . 

0 

is  the  state-bonus  modifier  matrix. 
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5 


■ 


I 


f = 


Fl.  KI2 


*12  K22 


m*m 


is  tiie  negative  of  the  positive  semi -clef  inite  cost  quadratic  form. 

The  matrix  entries  F..  are  functions  of  economic  parameters  and  the 

tj 

desired  terminal  inventory.  The  economic  significance  of  the  F..  are 

ij 

discussed  in  Munch  (1977),  pp.  18-19.  Refer  to  App.  B for  the  equa- 
tions defining  these  components. 


R = diag 


6,(1 


eif  11 ) ’ B5’ 


JQ  * 


is  the  bonus  cost  quadratic  form. 

F^  = lower  diagonal  m-1  * m-1  submatrix  of  F is  the  penalty  cost 
part  of  the  cost  quadratic  associated  with  the  reduced  state  L . 


0 


0 


13 


1 ,m 


0 


l + |™— *-]  *m-l 


is  the  penalty  cost  matrix  associated  with  B , 

* » I 

first  year  group  in  year  t.  F.  . , j = 2,  3 

entries  in  the  first  row  of  F. 


the  bonus  for  the 
m are  the  last  m-1 


DEFINITION  OF  THE  OPTIMAL  CONTROL  PROBLEM 

The  optimal  control  problem  is  to  find  a sequence  of  control  vec- 
tor (yearly  bonuses)  B such  that 

77  n 


) t (I  ,K  ) is  a minimum  subject  to: 

'■  n t t 


, = Al,  + 6B  + C , 
t + 1 I t 


4 - * 


a prescribed  terminal  inventory 


I - 1 / 2 B ' R B 
t t t 


1/2(L(  + SB ^ - L) ' + SB  - L) 


i B 1 1"  J ( I't  + S\  - <•>  • 


SOl.l'T ION  TO  TUT.  OPTIMAL  CONTROI.  PROBLEM 

Tbe  solution  to  the  optimal  control  problem  is  found  by  means 
of  the  discrete  version  of  the  maximum  principle.  Thus,  if  the  system 
Hamiltonian,  11^  = 11^  (l.f  , X | , B^  ) , is  defined  as: 


",  ■ -f,  * \’+i  (AI't  + + <:>  • 


then  the  optimal  solution  must,  at  each  stage  t,  satisfy 


3Ai+i  * Lt+1 


\ * 


(3) 


I 


Equation  (1)  is  seen  to  be  equivalent  to: 


-RBt  + S'F.,(I,  + SB(  - !.)  + B j F ' j ( I . t - I.) 


+ B . (F  S + S’Fj ) B + B'X. = 0 


Bt  = (R  - S'F2S  - F j ) ~ I (S  * F + BjF^I^  - 1)  + B'At+|]  , (4) 


wiie  re 


Fj  = B, (FjS  + S'Fj) 


Equation  (2)  will  yield  a re-statement  of  the  transition  equation 


constraint,  while  l.q.  (3)  yields: 


\ ■ >V'.t  + SB,  - I.)  + B.FjB,  * A'Xt+, 


Substituting  the  righthand  side  of  Eq . (4)  yields: 


A = F.,(L  - I.)  + A'A  .. 

t 2 t t+1 


+ (F2S  + BjFj)(R  - S'F2S  - F])“1|(S'F2  + B j F , ) ( L t - L)  + B'A^I 


Expanding,  collecting  terms,  and  solving  for  A j in  terms  of  1.^  and 


Vi  ■ ",  1-V-t  - '■>  * V • 


where 


F2  + (F2S  + BjF|)(R  - S'F2S 


F*>  ‘ (S’F2  + B , F t ) 
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and 


Bj  | A ' + (K,S  + B Fj )(R  - S'F,S 


B'  I • (')) 


Now  substituting  tin-  righthand  side  of  Eq . ( r> ) for  back 

into  Lei . (4)  and  substituting  tho  resulting  expression  lor  B in  terms 

of  a linear  combination  of  l.(  and  A(  , 


1.  . 

p ! 

Q 

i. 

t + 1 

i 

i 

i 

i 

t 

At+> 

i 

i 

i 

v ! 

i 

y 

A 

t 

- 

_ 

_ 

(10) 


where 


P = | A + (3(R  - S'  F.,S  - i'")  (S'F?  + B F - B'B,  A,)| 


1 I 


(ID 


Q = BI  (H  - S ' F ,S  - Fj‘)  ' 3'Bj 


(12) 


V - . 


(1  i) 


Z = BI'  • 


(14) 


I)  = 


* -1 


B(K  - S ' F.,S  - F t ) (S'F2  + BjFj  - B'B  ‘AjJL  + C 


(15) 


These  last  results  show  that  there  is  a discrete  linear  system 
whose  terminal  value  L , after  t steals,  depends  upon  the  value  of  the 
Initial  co-state  vector  A(j. 


-li- 


lt tlu’  a!  is  well  as  certain  other  model  parameters  are  nonzero, 
the  problem  is  well-posed  in  the  sense  that  all  the  matrices  that 
must  he  inverted  are  nonsingular. 

It  the  transition  matrix  of  Eq . (10)  is 


(16) 


then  it  is  easy  to  show  by  induction  on  n that 


- 

_ 

Lt+1 

Lo 

— 

A t+1 

= 

— 

A 

A„ 

t+1  J 

0 

+ (I  + <J>  + + ...  (f)^)!)  . (17) 


Assuming  that  if>  does  not  have  a unit  eigenvalue,  Eq . (17)  can  be  re- 
written as 


- 

“ 

Lt+1 

A t+1 

Lo 

= (P 

A 

A. 

t+1 

_ 

0 

+ (!-$)  *(1  - <t>t+1)D  . 


(18) 


In  addition  if  (}>  is  assumed  to  have  no  multiple  eigenvalues,  then 
4)  = V d lag  pjV"1  , (19) 


where  p . are  the  eigenvalues  of  (j>  and  V is  the  complex  matrix  whose 
columns  are  the  eigenvectors  of  <p.  Therefore,  Eq . (18)  can  be  written 
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■ 

- 

Ln+1 

,,  . ■ n+1  , 1 

'■<) 

=■  V diag  p.  V 

X 

X, 

t + 1 

() 

- 

- 

(20) 


+ V 


(1 


d i ag 


n+  1 
wi  ) 


(1  - M j ) 


V I) 


Tin-  de  i -iimpos  i t ion  of  P in  i'.tj . (19)  permits  the  potentially  numer- 

ically unstable  process  of  computing  the  powers  of  >p  to  be  replaced 
bv  the  more  stable  computation  of  powers  of  t lit'  eigenvalues  of  <p. 
Moreover,  till  spectral  decomposition  gives  some  insight  into  system 
stability  in  both  the  numerical  and  Lyapounov  senses.  Thus  the  sys- 
tem will  be  Lvapounov-stabl e and  numerically  stable  il  and  only  if 
tlu  eigenvalues  are  all  less  than  one  in  magnitude.  This  has  not  been 
tilt  case  with  any  of  the  models  tested  thus  tar.  Hence,  the  method- 
ology developed  here  will  fail  for  a sufficiently  large  problem  dimen- 


sion or  ext en t . 


If  the  terminal 

state 

and  co-state 

and  A.j . , then 

1 

' Lo  ' 

'T 

= 4> 

XT 

x<> 

a f t er  T 


+ D , 


years  are  denoted  as 


(21) 


where 


4i  = v 


diag 


v-'d 


(22) 


V 


-1 


D . 


i)  = V diag 


-I  1- 


I’art  i t inn  ing  f>  into  four  blocks  and  I)  into  upper  and  lower  parts, 


p-,2 1 r 'o 


V | ! *0 


't  ' K11L0  + I:12A0  + Du 


X0  = E12(LT  ~ E,1L0  - V 


provided  V.  | , is  nonsingular.  (The  nonsingularity  of  K ,}  is  a suffi- 
cient but  not  necessary  condition  for  the  solution  of  Eq . (24)  for  A(^.) 

Given  I.(J,  and  X as  determined  from  F.q . (25),  it  is  possible  to 
const  nut  each  state  and  co-state  from  the  preceding  according  to 


Lo  ■ 'o 


\)  ■ ’o 


“o  • W,|U2(L0  - I.,)  ♦ 


ho  • A,'l  + BB.  ♦ ■ 


Vi  ’ "i  |Aiat  - V *V  • 


Vi  ’ WVi  - 4*  + Vo,1  • 


(t  - I,  2 T-l ) , 


Having  sequentially  constructed  the  state,  co-state,  .and  optimal 
bonuses  as  Kqs.  (2b)  to  (28)  above,  it  is  possible  to  construct,  fo 
each  year  t,  an  augmented  and  modified  state  vector,  I.  , such  that 


m*  1 


KlBY< 


i.  + sit 

t n 


( 3 


The  rows  of  the  modified  (1-table  are  given  bv 


(1 

t , 1 


j 5 m 

j = m+1 , . . . , 2m-  1 


(3 


The  bonus  cost  , 


B , t 


of 


the  controlled  system  in  year 


is  given  bv 


(' 


( 3 


The  penalty  cost,  C of  the  controlled  system  is  given  by 
' % i 


= -1/2(1.  - 1,  ) ’ F’(L 


LJ  . 


Tot  .1 1 system  cost,  C.^.  { is  given  hy 


‘‘TOT,  t ('B,t  + l'P,t 


It  should  be  noted  that 


f<LfV“t>  ■ Sot., 


The  system  Hamiltonian,  li^  , in  year  t is  given  by 


Ht  ‘tot , t + Xt+iLt+i 


In  contrast  to  the  optimally  controlled  system,  the  uncontrolled 
system  evolves  according  to 


where 


L*+]  = A*L*  + IKj,  0 OP  , 


0 a. 


0 


0 


-lfl- 


Ml.  PROGRAM  OPERATING  INSTRUCTIONS 


FUNCTION 

The  DBMMI’  (Dvnamic  Bonus  Management  Model  Program)  computes  the 
solution  to  .1  discrete,  linear  optimal  control  problem  by  using  the 
discrete  maximum  principle.  The  problem  is  formulated  in  terms  of 
matrix  and  vector  algebra  and  a spectral  decomposition  of  the  basic 
state/co-state  transition  matrix  plays  a key  role  in  the  solution 
methodology.  The  cost  function  quadratic  form  matrix  entries  are  pre- 
computed as  functions  of  a set  of  economic  model  parameters. 

ARCHITECTURE 

The  DBMMI’  consists  of  a main  program  and  some  22  subroutines. 

All  are  coded  in  FORTRAN  IV  double  precision  for  the  IBM  370-158. 

Both  real  and  complex  arithmetic  are  used.  Total  source  program  length 
is  about  1700  instructions. 

ST(  >RA( ; i.  RE({U  I K KM KNTS  AND  TIM  INC  I 

Core  storage  requirements  for  the  program  and  its  subroutines  are 
dependent  on  problem  dimensions.  A sample  problem  characterized  by 
eight  year  groups  over  a ten  year  manning  horizon  required  22bK  bytes 
of  core  at  the  GO  step  with  object  code  produced  by  the  FORTRAN  G 
compiler.  This  problem  was  solved  in  a running  time  of  5.4  170-158 
CPU  seconds. 

LIMITATIONS 

Because  of  the  problem's  characteristic  stability  structure,  its 
dimensions  cannot  be  pushed  too  far.  Numerical  difficulties  were 
encountered  in  the  sample  problem  as  the  manning  horizon  was  varied 
significantly  beyond  ten  years. 

INPUT 

Table  2 lists  DBMMI’  Input  parameter  program  names  and  definitions. 
The  first  two  input  parameters  define  the  dimensions  of  the  optimal 
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Table  2 


DBMMP  INPUT 


Paramet  er 


Name 

[)f  f in  i t ion 

Uni  ts 

Or igin 

Nl 

Number  of  year 
groups  less  one 

Dimens i on  1 ess 
( i nteger ) 

Card 

i nput 

M 

Problem  extent 

Years  (integer) 

Card 

input 

T 

Problem  extent 

Years  (real) 

Card 

i nput 

IXZE.RO 

Initial  year  group 
state  vector  at 

t ime  zero 

Men  (real) 

Card 

input 

EI.T 

Terminal  year  group 
state  vector  aug- 
mented by  one  in 
dimension  to  include 
the  first  year  group 

Men  (real) 

Card 

input 

1)1 

First  distribution 
paramet er 

Dimension  1 ess 
( rea  1 ) 

Card 

input 

1)2 

Second  distribution 
parameter 

Dimens  ion  1 ess 
(rea  1 ) 

Card 

input 

W2 

Second  term  wages 

Dollars  ( rea 1 ) 

Ca  rd 

input 

Kilt) 

Subst i t ut ion 
paramet er 

Dimensionless 
( rea 1 ) 

Card 

input 

F.TA 

Demand  elasticity 

Dimension  1 ess 
(real) 

Card 

input 

EN 

Returns  to  scale 
parameter 

Dimensionless 

(real) 

Card 

input 

HP 

Supply  parameter 
vector 

Dimensionless 

(real) 

Card 

input 

AA 

No-bonus  accession 

rate 

Dimension  1 ess 
(real) 

Card 

input 

ALPHA 

No-bonus  continua- 
tion rate 

Dimension  less 
(real) 

Card 

input 

control  problem.  N1  is  tin-  number  ol  different  year  groups.  Currently 
Nl  must  always  be  set  to  7.  M,  the  problem  extent  or  manning  horizon, 
should  not  go  much  beyond  12  if  numerical  difficulties  are  It)  be 
avoided.  T is  just  a real  variable  convenience.  El, ZERO  is  the  Nl 
dimensional  vector  whose  components  are  the  number  of  men  in  year 
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groups  2,  3,  ....  Nl+1.  KLT  is  the  Nl  + 1 dimensional  vector  of  desired 
manpower  level  in  each  of  the  year  groups  1,  2,  ....  Nl+1.  The  aug- 
mented lirst  component  of  I,  is  the  constant  input  manpower  level  in 
year  group  1.  1)1,  the  distribution  parameter  (6.);  1)2  (6,J;  and  W2 

are  defined  in  Table  1.  RHO  is  p,  UTA  is  t,  UN  is  p,  HP  is  the  vector 
ot  8’,  AA  is  the  no-honus  accession  rate  equivalent  to  A,,  and  ALPHA 

1 I 

corresponds  to  u,  all  as  defined  in  Table  1. 

Table  3 defines  the  input  data,  card  formats,  and  deck  structure. 
Table  4 defines  the  dob  Control  Language  (JCI.)  card  requirements  for 
executing  a program  load  module  whose  Data  Set  Name  is  l.ATKST  belong- 
ing to  the  partitioned  data  set  PANBONUS. 

OUTPUT 

Because  the  current  version  of  the  DBMMP  is  a prototype  research 
tool,  much  of  its  output  is  generated  for  purely  diagnostic  purposes. 
Ibis  section  will  describe  only  the  most  important  output  parameter 
va  1 ues . 

Table  5 shows  a sample  output  echo  of  each  of  the  basic  input 
parameters.  N1  is  the  reduced  basic  problem  dimension,  N is  the  row 
and  column  dimension  that  will  be  dealt  with  in  the  solution  of  the 
control  problem,  M is  the  problem  extent  in  years,  and  T is  the  float- 
ing point  problem  extent. 

The  real  N1  * 1 matrix  LZRO  is  the  initial  vector  of  manpower 
level  in  year  groups  2-8.  The  real  N * 1 matrix  LT  is  the  desired 
final  vector  of  manpower  levels  In  year  groups  1-8.  The  bottom  row 
of  the  table  is  an  echo  of  the  economic  model  parameters  previously 
described . 

Table  6 is  the  N x 11  negation  semi-definite  matrix  associated 
with  the  penalty  cost  for  excesses  or  deficiencies  in  targeted  manpower 
levels.  It  is  computed  as  a function  of  LZF.RO,  LT,  and  the  parameters 
L)1  , 1)2  , etc  . 

Table  7 is  the  N * N state/co-state  transition  matrix  that  moves 
the  year  groups  and  associated  co-state  of  shadow  prices  from  one  year 
to  the  next.  It  is  Identical  to  the  matrix  defined  in  Kq . (16). 


I NPl'T  CARD  SPECIFICATIONS 
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Table  4 
PANBONUS  -1(1. 


//S4610P01  -IOB  ( 1191,  f»() ,51),  "FEBRUARY"  , CI.ASS=F 

//GO  EXEC  PGM* LATEST, REG  ION-2 30K 

//STEPLIll  1)1)  I)SN=S . S46  1 0 . A3  1 9 I . PANBONUS  ,1)1  SP=SHR 

/ /GO.  FT05F001  i)I)  I)DNAME=SYSI  N 

/ / GO  . FT06F00  1 1)1)  SYSOUT=A 

/ /GO . FT09F001  DD  SYSOUT=A,  1)CB=  (RECFM=FA,  I.RF.CL*  I 1 1 , BI.KS  I 7.E-  1 13) 
//GO. FT10F001  DU  SYSOUT=A,  DCB=  (RECFM=FA,  LRECL=  1 1 3 , BI.KS  IZE=  1 33) 
/ 'GO.  FT  1 1F001  1)1)  SYSOUT=A,  1)GB=  (RECFM*FA,  I.RECL*  1 3 3 , BI.KS  IZE=  1 33) 
/ /GO. FT12F001  DD  SYSOUT*>A , DGB* ( RECFM* FA , l.RECI.- 1 33.BLKSIZE-1 3 3) 
//GO.FT1  3F001  DI)  SYSOUT=A,I)CB=(RECFM=FA, I.RECL*  1 33,BLKSIZE=1  33) 
/ /GO . SYS  I N 1)1)  * 

7 H 8.0 

350.  250.  250.  160.  170.  180.  220. 

323.  307.  292.  278.  200.  200.  200.  200. 

.45  .55  10000.  .25  .68  .032  .005  323. 

.95  .95  .95  .72  1.0  1.0  1.0  .95 

/* 

/* 


Table  8,  (M+l)  x (2N-1),  Is  a tableau  of  the  year-by-year  optimal 
state/co-state  solutions.  The  lefthand  row  numbers  should  be  reduced 
by  one  to  give  the  appropriate  year.  The  first  eight  entries  in  each 
row  are  the  manpower  levels  in  each  year  group  i,  i = 1,  2,  ....  8, 
and  the  last  seven  entries  are  the  corresponding  co-state  variables. 

Table  9 summarizes  the  costs  associated  with  optimal  system  con- 
trol and  compares  them  with  the  costs  incurred  without  control.  To  be 
precise,  however,  the  year  In  the  second  column  at  the  left  should  be 
reduced  by  one.  PENALTY  COST  is  the  quadratic  component  of  the  cost 
associated  with  year  group  deviations  from  desired  final  manpower  levels. 
BONUS  COST  is  the  quadratic  cost  component  associated  with  bonus  pay- 
ments. HAMILTONIAN  is  the  discrete  system  Hamiltonian.  BONUS  1 and 
BONUS  5 are  the  optimal  bonuses  paid  to  men  entering  the  first  and 
fifth  year  groups.  PENALTY  COST*  is  the  pure  quadratic  penalty  cost 
corresponding  to  a no-bonus  policy.  NET  GAIN  is  the  difference  in 
the  cost  of  the  uncontrolled  and  optimally  controlled  systems.  The 
last  line  of  the  table  totals  the  various  cost  components. 


DBMMP  INPUT  PARAMETER  ECHO 
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IV.  I’ ROM  RAM  DFSOKI  I’TION 


I’ROOKAM  FUNCTION 

Hu  Dvn.imit  Bonus  Management  Model  I'rogram  is  designed  t < • deter- 
mint'  .in  optimal  polity  tor  the  payment  of  yearly  bonuses  to  military 
personnel  in  a number  of  different  "year  groups."  The  bonuses  are 
treated  as  the  rontrol  variables  in  a discrete  linear  optimal  control 
system  in  which  the  components  of  the  current  state  vector  are  the 
number  of  men  in  each  of  a specif ied  number  year  group.  The  objective 
of  tlie  optimal  control  policy  is  to  drive  the  system  from  a given 
initial  state  in  year  zero  to  a prescribed  terminal  state  over  a spec- 
ified number  of  years,  while  minimizing  a given  cost  function.  This 
cost  function  is  the  summation  over  each  year  of  a positive  son i - 
definite  quadratic  form  in  both  the  control  variables  and  the  differ- 
ence between  the  year's  state  and  the  desired  terminal  state.  The 
quadratic  form  itself  is  computed  as  a function  of  a set  of  system- 
characterizing  economic  model  parameters  and  the  terminal  state. 

The  solution  for  the  optimal  policy  is  derived  as  a consequence 
of  the  disc  rete  version  of  the  maximum  principle.  essentially,  it 
reduces  to  the  construction  of  a discrete  linear  system  of  twice  the 
dimension  of  the  year-group  state  and  the  subsequent  determination  of 
an  unknown  initial  co-state  vector  that  will  yield  the  prescribed 
terminal  state. 

• In  addition  to  computing  the  sequence  of  yearly  bonuses  and  asso- 
ciated costs  for  the  optimal  bonus  policy,  the  program  also  performs 
a spectral  analysis  of  the  system  state/co-state  transition  matrix. 

It  also  produces  a net-gain  comparison  between  the  optimal  policy  and 
the  no-control  policy  of  zero  bonuses. 

I’ROCKAM  ARCH  I T MCI  l ' Rif 

The  DBMM1’  is  organized  into  a FORTRAN  IV  main  program  and  some 
22  subroutines.  Communication  between  main  and  subroutine  and  from 
subrout ine-to-subrout ine  is  accomplished  by  the  use  of  an  extensive 
common  block  structure  and,  where  appropriate,  calling  sequence 
parameters. 
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The  main  program  is  responsible  for  calling  the  data  input  and 
problem  setup  subroutine,  PROGO,  solving  for  the  required  initial  co- 
state  vector,  sequentially  computing  the  optimal  bonuses,  and  calling 
the  computational  and  analysis  subroutine  COSTS  and  OUTPUT  to  print 
out  the  states,  co-states,  bonuses,  associated  costs,  system  Hamil- 
tonians, uncontrolled  system  costs,  and  net  gains. 

The  subroutines  are  responsible  for  1/0  and  the  linear  algebra 
and  eigenanalysis  functions  that  solve  the  optimal  control  problem. 

PROGRAM  INPUT 

Program  input  is  handled  entirely  by  the  subroutine  PROGO.  For 
detailed  input  data  card  specifications  see  Sec.  III. 

PROGRAM  OUTPUT 

Program  outputting  functions  are  shared  among  the  main  program 
and  a number  of  major  subroutines.  In  its  current  research  tool  state, 
much  of  the  program's  output  is  diagnostic  and  the  reader  need  not  be 
concerned  with  it.  The  important  outputs  are:  the  echo  of  the  input 
parameters;  the  negative,  semi-definite  cost  quadratic  form,  F;  the 
system  state/co-state  transition  matrix,  R;  the  array  of  modified  year 
group  levels  (G  TABLE);  and  the  final  array  of  penalty  cost,  bonus  cost, 
total  cost,  system  Hamiltonian,  yearly  bonuses,  penalty  cost  of  the 
uncontrolled  system,  net  gain,  and  cumulative  cost  and  net  gain. 

PROGO  echoes  out  input  parameters  and  prints  out  the  matrix  of 
the  cost  quadratic  form.  The  PROGO-called  subroutine,  COMPUT,  writes 
out  the  system  state/co-state  transition  matrix.  The  main  program 
prints  out  the  transition  matrix  eigenvalues.  The  subroutine  COSTS 
writes  out  the  two  parts  of  the  A-Matrix,  and  OUTPUT  is  responsible 
for  the  unmodified  L and  ELAMB  table  output.  COSTS  also  handles  the 
outputting  of  cost/benefit  performance  parameters  and  the  yearly 
Hamiltonian. 

For  detailed  definitions  of  the  output  parameters,  see  the  ex- 
ample of  output  given  in  Sec.  Ill  and  the  description  of  the  relevant 
outputting  subroutine. 
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PROG  RAM  (’.( )M PUTATK )N  A I . S , K(£ll IE , NC  E 

Program  operation  begins  with  setting  the  maximum  of  the  size  of 
the  problem  dimension  to  AO.  PROGO  is  then  called  to  read  input 
parameters  and  set  up  the  various  matrices  and  vectors  that  character- 
ize the  problem.  PROGO  calls  COMPUT  to  generate  the  uncontrolled 
system  state  transition  matrix  and  cost  quadratic  form  matrix.  From 
these,  together  with  the  control  variable  matrix  multiplier,  PROGO 
generates  the  state/co-state  transition  matrix.  It  also  sets  the  com- 
plex part  of  the  transition  matrix  to  zero.  When  PROGO  returns  control 
to  the  main  program,  various  loop  limits  are  set  and  the  real  part  of 
the  state/co-state  transition  matrix  is  saved.  Next,  CBAL  is  called 
to  balance  the  state/co-state  transition  matrix,  COMHES  is  called  to 
compute  the  upper-Hessenburg  matrix  similar  to  it,  C0MLR2  finds  its 
complex  eigenvalues,  and  CBABK2  back-transforms  the  corresponding 
complex  eigenvectors  to  those  of  the  original  state/co-state  matrix. 

The  complex  matrix  of  these  eigenvectors  is  used  subsequently  to  sim- 
plify and  stabilize  the  required  computation  of  the  partial  sum  (up 
to  "T"  years)  of  a geometric  power  series  in  the  state/co-state  transi- 
tion matrix. 

The  resultant  geometric  series  partial  sum  matrix  is  then  parti- 
tioned along  with  initial  and  desired  final  states  to  determine  the 
required  initial  value  of  the  co-state  vector.  Once  this  initial  value 
is  determined,  the  main  program  computes  each  year  group  state  vector, 
co-state  vector,  and  optimal  bonuses  recursively  from  the  previous 
year's  values.  This  process  also  permits  the  sequential  construction 
of  each  of  the  output  modified  state,  co-state,  and  cost  analysis 
tabular  entries.  The  modified  state  table  is  used  both  as  a user  out- 
put and  computational  convenience. 

PROGRAM  FLOW  DIAGRAM  AND  SUBROUTINE  CALLING  MAP 

To  supplement  this  prose  description  of  the  program  computational 
sequence.  Fig.  1 presents  a flow  diagram  of  the  DBMMP  main  program. 

For  a definition  of  the  mathematical  symbols  and  equations  development , 
refer  to  Sec.  II. 
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l--Flow  diagram  of 


thf  DBMMP  main 


program 


Fig.  1 — continued 
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As  a furtl  r aid  to  understanding  tlie  program's  functional  flow, 
Kig.  2 gives  a si. 'routine  calling  map,  which  is  a compact  representa- 
tion of  the  program's  functional  flow.  Every  subroutine  appears  in 
order  of  its  first  use  by  the  calling  main  program  or  subroutine.  But 
the  subroutines  do  not  appear  again  if  they  are  called  repeatedly  by 
the  calling  subroutine  or  main  program. 

Si: B Kl JUTINK  .1  >ESCK  1 PT ION 

In  the  sequel,  each  subroutine's  full  input  is  described  in  terms 
of  explicit  calling  sequence  structure  and  those  parameters  that  are 
implicitly  accessed  through  COMMON. 

Subroutine  J’KOGO  (T,M, JK) 

PROGO' s primary  functions  a.e  total  program  input  and  construction 
of  problem-characterizing  matrices  and  vectors.  It  reads  in  the  dimen- 
sion of  the  problem  (number  of  different  year  groups  beyond  the  first 
year  group),  the  problems  extent  (number  of  transitional  years),  the 
initial  state  vector,  the  desired  terminal  state  vector,  and  the  prob- 
lem characteristic  economic  parameters.  From  these  are  computed  the 
entries  of  the  cost  function  quadratic  form.  In  addition,  PROGO  ini- 
tializes to  zero  the  real  and  imaginary  parts  of  the  state/co-state 
transition  matrix  whose  dimension  is  twice  the  dimension  of  the  number 
of  year  groups.  The  entries  of  the  matrix  are  filled  in  by  the  PROGO- 
called  subroutine  COMPUT. 

Refer  to  Table  3 for  PROGO's  input  parameter  names  and  definitions. 
All  PROW)  input  originates  by  card  input. 

Table  10  lists  the  set  of  PROGO  output  parameters,  definitions, 
and  units. 

On  initiation,  PROGO  sets  the  maximum  allowable  state/co-state 
dimension  to  40.  Next  it  reads  Nl,  the  basic  problem  dimension;  M, 
the  integer  problem  extent;  and  T,  the  real  floating  point  problem 
extent.  It  then  sets  to  zero  intermed iat e storage  and  the  real  and 
imaginary  parts  of  the  sta  /co-state  transition  matrix.  Then,  the 
initial  and  desired  final  state  vectors  are  read  in.  The  economic 
model  parameters  are  read  in  along  with  the  year  group  transition 
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MAIN 
PRO  GO 
PRINTR 
COMPUT 
PRINTR 
MAT MAT 
MATR1N 
TRANSP 
MATVEC 
CBAL 
COMHES 
COMLR2 
CBABK2 
PRINTC 
MAT RID 
MATDIA 
MATMTC 
MATVEC 
MATRIN 
GFIX 
MATMAT 
MATMAT 
OUTPUT 
COSTS 
GFIX 
STAR 

EXIT  (System  Subroutine) 


Fig.  2 — DBMMP  calling  map 
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Table  10 
PR0G0  OUTPUT 


Pa  rami' 1 < r 

Del  in i t ion 

Units 

r 

Negative  semi -def in i t e cost 
quad ra lie  f orin 

2 

Dollars/men  (real) 

A 

State  transition  matrix 

Dimens  inn  1 ess  (rea 1 ) 

BB 

Control  variable  matrix 
mu  1 1 ip  1 ier 

Men/dollars  (real) 

B 

Constant  vector  in  state/ 
co-state  transition 
equat i on 

(Men,  men/dollars | 
(real) 

AR 

Real  part  of  the  state/ 
co-state  transition 
ina  t r i x 

(a) 

A I 

Imaginary  part  of  the  state/ 
co-state  transition  matrix 

(a) 

BTZ 

Auxiliary  matrix  for  compu- 
tation of  stagewise  bonuses 

Men/ doll ars  ( rea  1 ) 

aBoth  AR  and  AI  have  the  units  structure 


Dimensionless  ! Dollars 

i 


Dollars  i Dimensionless 

i 


parameters.  Subsequently,  the  cost  quadratic  form  matrix  entries  are 
computed  as  a function  of  the  economic  model  parameters  1)1,  1)2,  W2, 
....  etc.  Once  F has  been  computed,  PROGO  branches  to  COMPUT  for  the 
computation  of  the  last  six  of  PKOGO' s output  parameters.  A,  BB,  B, 

AK,  AI,  BTZ. 

S u hr  out i ne  I ’ R I NT  R ( ( )UTMAT,N,M,MA , K ) 

PRINTR  prints  out  a real  matrix  according  to  a number  of  variable 
formats.  The  particular  format  is  selected  on  the  basis  of  an  input 
parameter.  Table  11  lists  PRINTR  input  parameters,  definitions,  and 
origins. 
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Table  1! 
I’HINTK  INPU'I 


Paramet er 

Del  ini t ion 

U nits 

Origin 

OUTPUT 

Location  of  the  first  entry 
of  the  matrix  to  be  printed 

Dimensionless 
( i nteger ) 

Ca  ! I i ng 
program 

N 

Row  dimension  of  the  matrix 

Dimension  1 ess 
(integer) 

Cal  ling 
program 

M 

Column  dimension  of  the 
m.i t r i x 

Dimension  less 
(integer) 

Ca 1 1 i ng 
program 

MA 

Maximum  row  dimension 

Dimens  ion  less 
(integer) 

Ca l 1 ing 
program 

K 

Heading  spei  ification 

Dimens i on  less 
( integer) 

Ca 1 1 ing 
program 

PRINTR's  output  is  row-by-row  entries  of  the  input  ra.it rix,  OUTMAT, 
ten  row  entries  of  1)12.4  format  per  line. 

Subroutine  COMPUT 

COMPUT  sets  up  the  basic  matrices  involved  in  the  optimization 
problem.  Table  12  lists  COMPUT's  input  parameters,  definitions,  units, 
and  origins.  Table  13  lists  its  output  parameters,  definitions,  and 
un its. 

COMPUT  begins  its  operation  by  setting  the  maximum  dimension  to 
Nl,  the  basic  overall  problem  dimension.  It  then  sets  to  zero  the 
effective  constant  input  vector,  CC,  and  the  state  transition  matrix, 

A,  and  several  auxiliary  matrices  S,  BB,  ST,  and  PI.  Nonzero  entries 
are  next  computed  for  the  state  transition  matrix.  A,  the  control 
variable  matrix,  BB,  and  the  effective  reduced  dimension  constant 
input,  CC.  The  positive  definite  quadratic  form  in  the  control  vari- 
ables, RR,  is  then  set  to  its  initial  value.  Subsequently,  the  aux- 
iliary matrix,  Fj,  is  constructed;  and  the  reduced  order  cost  quadratic 
format,  F^,  is  set  up  as  the  lower  Nl  * Nl  diagonal  block  of  the  orig- 
inal quadratic  form.  In  the  final  sequence  of  computations,  each  of  the 
four  Nl  x Nl  blocks  of  the  real  part  of  the  state/co-state  transition 
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Tab  It*  12 
COM PUT  INPUT 


Pa  rame  t t*r 

Do  f ini t i on 

Un  i t s 

Origin 

N 1 

Problem  dimension 

Dimension  1 ess 
( integer) 

PROGO 

ALPHA 

State  transition  parameters 

Dimens i on  1 ess 
( integer) 

PROGO 

HP 

Supply  parameter  vector 

Do  1 1 a rs~ 1 
( rea 1 ) 

PROGO 

AA 

Dimens i on  1 ess 
(rea  1 ) 

pr<  x ;o 

F 

Cost  quadratic  form 

Do  1 1 ars/  men 
(rea  1 ) 

PROGO 

Table  1) 

COMPUT  OUTPUT 

Pa r ame  ter 

Def ini t ion 

Un  i t s 

AR 

The  real  part  of  the  state/ 
co-state  transition  matrix 

(See  Sec.  1 1 ) 

A I 

The  zerti  imaginary  part  of 
the  state/co-state  transi- 
tion matrix 

(See  Sec  . 11) 

B 

The  constant  state/co-state 
input  vector 

(Men,  men/dol lar) 

BB 

The  control  variables  multi- 
plier matrix 

Men/do  liar 

BTZ 

Auxiliary  matrix  for  sequen- 
tial bonus  computation 

Men*" /do liar  (real) 

- )H- 


matrix,  AK,  is  computed  along  willi  the  constant  state/co-state  transi- 
tion-input vector,  it,  and  tin*  auxiliary  bonus  computation  matrix  BT7-. 

For  a detailed  understanding  of  the  meaning  of  PROOO's  operation  and 
symbols  see  Appendix  A. 

MATMAT,  the  DBMMP  System  Real  Matrix  Multiplication  Subroutine 

Table  1A  shows  MATMAT's  input  parameters,  definitions,  and  origins. 
A " " in  tin1  units  column  means  the  units  are  determined  by  the 

particular  input. 


Table  1A 
MATMAT  INPUT 


Parameter 

Def in i t ion 

Units 

Origin 

AMATIN 

Pro-multiplier  matrix 

(real) 

Ca 1 1 i n g 
program 

BMATIN 

Post  -mu  1 1 i p 1 ie  r matrix 

( rea  1 ) 

Ca  Ming 
program 

NI 

Row  dimension  of  /AMATIN 

Dimens i on  less 
( integer) 

Cal  ling 
program 

NK 

Column  dimension  of 

AMATIN  = row  dimension 
of  BMATIN 

Dimension  less 
(integer) 

Ca fling 
program 

N.l 

Column  dimension  of 

BMATIN 

Dimensionless 
( integer) 

Ca 1 1 i ng 
program 

The  output  of  MATMAT,  OUTMAT,  is  the  N1  * N.l  dimensional  product 
AMATIN  * BMATIN . 

On  entrance,  the  product  matrix,  OUTMAT,  Is  set  to  zero.  Then 
the  product,  according  to  the  usual  matrix  mul t Ipl icat ion  rules,  is 


cik ' • 


m being  the  common  column  dimension  of  (a^)  and  the  row  dimension  of 


(bik>- 
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Subrout  ine  MATKIN,  the  DBMMP  System  Real  Matrix  Inversion  Subroutine 
MATKIN's  input  parameters,  definitions,  units,  and  origins  are 
g i veil  i n I able  IS. 


Table  IS 
MATKIN  INPUT 


Par a me t e r 

Def in i t i on 

Un  i l s 

Origin 

A 

Matrix  to  be  inverted 

( rea  1 ) 

Cn 1 l i ng 
program 

N 1 

Maximum  row  dimension 

of  A 

( Tea  1 ) 

Ca  Iling 
program 

N 

Actual  row  dimension 
of  A 

(real) 

Ca  I M ng 
program 

B 

Matrix  of  right  band 
side  vectors  in  the 
problem  AX=B 

(real  ) 

Cn  Ming 
program 

Ml 

Kow  dimension  of  K 

Dimensionless 
( int eger ) 

Cal  1 i ng 
program 

M 

Column  dimension  of  B 

1)  i mens  ion  1 ess 
( integer  ) 

Ca Iling 
program 

I PIVOT 

Internally  used  one- 
dimension  integer 
array 

Dimensionless 
( i nteger ) 

Ca 1 1 i ng 
program 

INDEX 

Internally  used  two- 
dimensional  integer 
ar  ray 

Dimens  ion  1 ess 
( integer) 

Ca 1 1 ing 
program 

ISINO 

Singu 1 ar /nons ingu lar 
flag 

Dimensionless 

(integer) 

Ca  Ming 
program 

DKTEKM 

Determinant  of  A 

(real  ) 

Ca Iling 
program 

MATKIN's  output,  A * , is  stored  in  A.  A is  destroyed  and  must 
be  saved  before  the  call  to  MATKIN  if  it  is  to  be  used  again.  MATKIN 
uses  Ciauss-. Iordan  pivotal  elimination  with  full  pivoting.  At  each  of 
N I stages,  the  element  of  maximal  absolute  value  is  seleeted  from  among 
those  row  and  column  entries  that  have  not  previously  been  used  as  a 
pivot  element.  Rows  are  then  permuted  to  bring  the  pivotal  element 
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into  thi*  diagonal  position  and  the  column  entries  above  and  below 
the  new  maximal  element  are  annihilated  by  appropriate  elementary  tow 
operations.  The  row  interchanges  are  kept  track  of  and  at  the  end 
of  the  last  stage  the  columns  of  the  inverse  are  permuted  inverse  to 
the  original  row  transpositions  to  unscramble  them. 

TRANS P } the  DBMMP  System  Matrix  Transposition  Subroutine 

TRANSP's  input  parameters,  definitions,  units,  and  origins  .ire 
given  in  Table  16.  Its  output,  XOUT,  is  the  N.I  x N1  transpose  ot  XIN. 


Table  16 
TRANSP  INPUT 


Parameter 

Def in i t ion 

Units 

Origin 

XI N' 

The  matrix  to  be 

(real ) 

Cal  1 ing 

transposed 

program 

NI 

Row  dimension  of 

XIN 

Dimens  ion  1 ess 

Ca 1 1 ing 

( integer ) 

program 

N.I 

Column  dimension 

of  XIN 

Dimension  1 ess 

Call ing 

( integer) 

program 

MATVEC , the  DBMMP  System  Matrix,  Vector  Multiplication  Subroutine 

MATVEC's  input  parameters,  definitions,  units,  and  origins  are 
given  in  Table  17.  MATVEC's  output  is  the  N dimensional  real  vector 

Table  17 
MATVEC  INPUT 


Parameter 

Def in i t ion 

Uni  ts 

Origin 

AMATIN 

Pre-mu  1 1 i p 1 y i ng  matrix 

(real  ) 

Ca 1 1 ing 
program 

VF.CTIN 

Post-multiplying  vector 

(real  ) 

Ca 1 1 ing 
program 

Nl 

Maximum  row  dimension 

of  AMATIN 

Dimension  1 ess 
( integer) 

Ca  1 1 i ng 
program 

N 

Row  dimension  of  A 

Dimensionless 

(integer) 

Ca 1 1 i ng 
program 

J 
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VLU'UT.  At  tor  sot  tint;  t lit'  product  vector,  VLCUUT , to  zero,  MATVLC 
■ omputcs  t ho  m.it  r i x-vec  t or  product  according  to  the  usual  rule: 


b. 

l 


l 

j=> 


a . . x . . 

• I I 


Subrout  inc  l BAl . (NM , N , AK  , A I , l.t  >_W  , 1 ( dl , SCAI.I! ) 

CBAL  is  used  to  pro-condition  the  complex  state/co-state  transi- 
tion matrix,  AK  + i A I , (AI=0),  balancing  it  and  isolating  its  eigen- 
values whenever  possible.  Balancing  is  performed  by  elementary  row 
and  column  interchanging  and  scaling.  CBAL  is  part  of  the  Arginine 
National  Laboratory  1. 1 S PACK  eigenvaluo/eigenvector  package.  Real 
arithmetii  is  used  throughout  t.'BAL. 

LBAl.’s  input  parameter  definition,  units,  and  origins  are  given 
in  Table  18.  CBAL's  output  parameter  del inition  and  units  are  shown 
in  Table  Id. 


Table  18 


CBAL  INPUT 


i’arameter 

l)e  f ini  t ion 

Uni  ts 

Origin 

NM 

Maximum  row  dimension  of 
matrix  to  be  balanced 

Dimension  1 ess 
(integer) 

Ma  i n 

N 

The  order  of  the  matrix 

Dimensionl ess 
(integer) 

Ma  i n 

AK 

Heal  part  of  the  matrix 

See  Sec . I I 

Ma  i n 

A I 

Imaginary  part  of  matrix 

See  Sec . 1 1 

Main 

Subrout  ine  ( X )M  1 i KS  (MA , N,AK  ,AI,1 .0 , IGli  ,1  NT) 

DBMMP  uses  COMHKS  to  reduce  the  balanced  matrix  output  of  CBAL 
to  upper  Hessenberg  form  by  stabilized  elementary  similarity  trans- 
formations. COMHKS  input  parameter  definitions,  units,  and  origins 
are  shown  in  Table  20. 
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Table  19 
CBAI.  OUTPUT 


P.l  I'.IIIU-  III 

I)i’  1 i n i t ion 

Un  i t s 

AH 

The  real  part  of  the 
balanced  matrix 

See  Set  . 1 1 

A I 

The  imaginary  part  of  t lie 
balanced  matrix 

See  Sec.  II 

1.0 

Indices  specifying  location 
of  zero  block  in  the 
balanced  matrix 

1)  i mens  ion  1 e: 
( integer ) 

IUH 

Indices  specifying  location 
of  zero  block  in  the 
balanced  matrix 

1)  i mens  ionic: 
( integer) 

SC.AI.F. 

N dimensional  vector  of 
scaling  information 

1)  i mens  ion  1 e: 
( rea 1 ) 

Table  20 
COMHKS  INPUT 


Parameter 

De f i n i t i on 

Un  i t 

Or  igin 

MA 

Matrix  row  dimension 

of  the  balanced  matrix 

Dimension  1 ess 
(Integer ) 

Ma  i n 

program 

N 

Kow  dimension  of  the 
balanced  matrix 

l)i mens  ion  1 ess 
( integer ) 

Ma  i n 

program 

LOW  1 

u;h  f 

Balanced  matrix  zero 
block  bounding  indices 

I) i mens  i on  less 
( integer) 

Ma  i n 
by 

program 

CBAI. 

AH 

Real  part  of  the 
balanced  matrix 

See  Sec . I I 

M a i n 
by 

program 

CBAI, 

A I 

Real  part  of  the 
balanced  matrix 

See  Sec . 1 I 

Ma  i n 
by 

program 

CBAI. 
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Out  put 


COMHKS' s output  parameter  definitions  and  units  are  given  in 


Tab  1 1 


Table  21 


COMHKS  OCTI’UT 


amel 

or  l)ef in i t ion 

I'n  i t s 

AK 

Rea)  part  of  the  upper 
liessenberg  matrix 

See  Sec . 

1 1 

A 1 

Imaginary  part  of  the 
upper  liessenberg  matrix 

See  Sec . 

1 I 

I NT 

Integer  array  containing 
row— columns  inter- 

Dimens i on 
( i nt  ege 

less 
r ) 

change  information 


Sul. rout  ine  C.t >MI  .R  2_  ( MA  , N ,1,0,1011,1  NT  , AK , A I , WR , W 1 , X K , Z I , 1 IKK ) 

C0MLR2  computes  the  eigenvalues  and  eigenvectors  of  t lie  reduced 
liessenberg  matrix  output  of  COMHKS.  OOMLR2 ' s input  parameter  del ini- 
t ions,  units,  and  origins  are  shown  in  Table  22.  C0!U,K2's  output 

parameter  definitions  and  units  are  given  in  Table  21. 

Subrout  i m CBABK2  (MA , N , 1.0 , I OH , SOAI.I. , N , ZR  , Z I ) 

CBABK2  back- 1 rails  forms  the  eigenvectors  >>1  tile  upper  liessenberg 
matrix  obtained  from  the  output  of  C0MI.R2  to  the  eigenvectors  of  the 
original  unbalanced  state/co-state  transition  matrix. 

0BABK2's  parameter  definition,  units,  and  origins  are  shown  in 
Table  24.  On  output  ZR  and  ZI  are  the  real  and  imaginary  parts 
of  the  eigenvectors  of  the  original,  complex,  unbalanced  state/co- 
state matrix. 

1’RINTC,  the  Com])  1 ex  Analog  of  PRINTR 

PRINTO  prints  out  a complex  matrix  preceded  by  an  accompanying 
variable  format.  Its  input  parameter  definitions,  units,  and  origins 
are  shown  in  Table  25.  PKINTC's  output  is  variable  header  101)12.4 
format  with  five  complex  entries  (real,  imaginary)  per  output  line. 
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T.il.lr  22 


COM1.K2  I NI’UT 


Paramot  el 

Del  in i t inn 

Un  i t s 

Origin 

MA 

Max i mum  d i mens  ion 
nf  the  Hessenberg 
ma t r i x 

1)  i mens  i on  1 ess 
( integer ) 

Ma  i n 
by 

program 

COMIIKS 

N 

Dimension  of  the 
Hessenberg  matrix 

I)  i mens  ion  1 ess 
( i nteger ) 

Ma  i n 
by 

program 

COMHKS 

MOW 

Zero  block  boundary 
indices  determined 
by  CBAI. 

IGH 

Zero  block  boundary 
indices  determined 
by  CBAI. 

I NT 

Vector  of  row-column 
permutat ions 

AK 

Real  part  of  the 
upper  Hessenberg 
matrix 

See  Sec . 1 I 

A 1 

Imaginary  part  of  the 
upper  Hessenberg 
matrix 

See  Sec  . I I 

Table  23 
C0MLR2  OUTPUT 


Parameter 

Def ini t ion 

Units 

WR 

Real  part  of  the  eigenvalues  of 
the  upper  Hessenberg  matrix 

Dimension  1 ess 
(real) 

WI 

Imaginary  part  of  the  eigenvalues 
of  the  upper  Hessenberg  matrix 

Dimension  1 ess 
(real) 

ZR 

Real  part  ot  the  eigenvector 
of  the  upper  Hessenberg  matrix 

(Men,  men /do liar] 
(real) 

ZI 

Imaginary  part  of  the  upper 

Hessenberg  matrix 

[Men,  men/dollar] 
(real) 

I ERR 

Error  i nd  1 ca tor : 

0 if  all  eigenvalues  determined, 
j It  jth  eigenvalue  not  determined 

Dimensionless 
( integer) 
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Tab  1 1 

• 24 

CBABK2 

INPUT 

Paramet  er 

l)ef  ini  t ion 

Un  i t s 

Origin 

MA 

Maximum  dimension  of 
original  unbalanced 
ma t r i x 

Dimens  ion  1 ess 
( int eger ) 

Ma  i n 

program 

N 

Dimension  of  the  matrix 

Ma  i n 

program 

LOW  ) 

Kill  f 

Zero  block  boundary 
i nd  i i'i's 

Ma  i n 
by 

program 

CBAL 

scam: 

Row-column  formula! ion/ 
scaling  information 

(rea 1 ) 

Main 

by 

program 

CBAL 

ZR 

Real  part  of  the  eigen- 
vector of  the  reduced 
llessenberg  matrix 

(Men,  men/dollar] 

Main 

by 

program 

COMLRZ 

Zl 

Imaginary  part  of  the 
e i genvec tor 

(Men,  men/dollar] 

Ma  i n 
by 

program 

COMLKZ 

Tab  1 e 2 4 


PH  INTO 

INPUT 

I’aramet  er 

Def in i t ion 

Units 

Or  igi  n 

OUTMAT 

Complex  matrix  to 
to  be  printed  out 

(comp  1 ex ) 

Ma  i n 

program 

N 

Row  dimension  of 

OUTMAT 

(real) 

Main 

program 

M 

Column  dimension 
of  OUTMAT 

(real) 

Main 

program 

MA 

Maximum  row 

( rea  1 ) 

Ma  i n 

program 

dimension  of 
OUTMAT 
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MATKID,  DBMMP's  Complex  Matrix  Inversion  Subrout  ine 
Analog  ol  MAI K IN 

Set ■ subrout ine  MATRIN  for  suhront  ine  railing  sequence,  I/O  param- 
eter ile  t i n i t i ons  , . omput  at  ion.il  sequence,  etr. 

Suhrout ine  MAID  I A (AMAT I N , D I AC  I N , N , N 1 ,OUTMAT) 

MATDIA  pos t -nut  I t i p 1 i es  a complex  square  matrix  by  a complex  diag- 
onal matrix.  MATDIA' s input  parameter  definitions,  units,  and  origins 
are  shown  in  Table  2b. 


Table  2b 
MATDIA  INPUT 


Parame t e r 

Def in i t ion 

Tn  i t s 

Or  i g i n 

AMAT 1 N 

Pre-multiplier  matrix 

( comp  lex) 

Ma  i n 

program 

D I AC  IN 

Vector  of  post- 
mu  1 t i [>  1 i cr  m.i  t r i x 
diagonal  entries 

( comp  1 ex ) 

Ma  i n 

program 

N 

Dimension  of  AMATIN 

Dimens  ion  1 ess 
( integer ) 

Ma  i n 

program 

N1 

Maximum  row  dimension 

Dimens i on  1 ess 
( integer) 

Ma  i n 

program 

* 

MATDIA  s output  is  the  complex  product  AMATIN  diag  (DIACIN(.I) ) . T.ach 
entry  in  the  jth  column  of  AMATIN  is  multiplied  by  DIACIN(.J)  . 

MATMTC,  the  DBMMP  Systems  Complex  Matrix  Multiplication  Subroutine 
MATMTC  functions  analogously  to  MATMAT.  For  calling  sequence 
parameter  definitions,  etc.,  see  section  MATMAT. 

Subroutine  CFIX  (M.ITIMK) 

CFIX  modifies  the  state/co-state  versus  year  table  so  that  it  can 
be  used  by  the  costs  subroutine.  It  computes  and  inserts  the  leading 
year  group  and  entry  and  modifies  the  fifth  year  group  entry.  In 
addition,  CFIX  computes  the  matrix  constant  used  in  the  sequential 
computation  of  the  optimum  bonuses.  CFIX  input  parameter  definitions. 


-47- 


units,  ci  nd  origins  are  shown  in  Table  27.  OK  I X output  parameter 
<ie!  init  ions  and  units  are  shown  in  Table  28. 

Tali  I e 27 
(■FI  X INPUT 


i rand  e i 

Def init i on 

Unit  s 

Origin 

M 

Problem  extent 

Years  (ini eger ) 

Main  program 

1 T I Ml 

Year  flag 

1)  i mens  ion  (ini  eger  ) 

Main  program 

N 

l ot  a 1 prob  1 eni 
d i mens  ion 

Dimens  ion  ( int  eger ) 

Main  program 

BTZ 

Aux  ilia  rv  mu t r i x 

> 

Men  /do  liar  ( rea 1 ) 

(lOMPl'T 

A 1 

Auxi 1 i ary  mat r i x 

Dimension  (real) 

OOMPl'T 

PI 

Auxi 1 iary  mat r ix 

2 

Do  1 1 ars/ men  ( rea 1 ) 

COMIT'I 

HP 

Supply  parameter 
vec  t or 

Do  liars  * ( rea  1 ) 

I’KOOO 

IT 

Bonus  values 

Do] 1 ars  ( rea  1 ) 

Main  program 

On  entrance  ('.FIX  checks  the  ITIMK  input  flag  to  determine  whether  it 
is  the  final  time  through.  If  ITIMK  = 1 and  it  is  the  first  time 
through,  OF  I X computes  the  looping  indices  NPt.USl  and  MPI.TS1  and  the 
auxiliary  matrices  BTZA  and  TF.RM1  , which  will  be  used  by  the  main  pro- 
gram for  sequential  bonus  computation  and  returns  to  the  calling  pro- 
gram MAIN.  Otherwise,  looping  over  each  year,  it  modifies  the  first 
and  fifth  entries  of  the  ('.-table  for  the  given  year  specified  by  the 
value  of  M by  adding  linear  multiples  of  the  optimal  bonuses  for  that 

Table  28 
OF  IX  OUTPUT 

Parameter  Definition  Units 

BTZA  Auxiliary  matrix  Men* /dollar  (real) 

TF.KM1  Auxiliary  matrix  Dimensionless  (real) 

0 Modified  O-table  |Men,  men/dollar|  (real) 
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vt'.ir.  At  the  com- 1 iis  ion  ol  the  loop,  CFIX  returns  to  the  cal  ling 
program,  which,  in  this  case,  is  the  subroutine  COSTS. 

Subrout  in.  OF  1 ITT  (.IK  ,.l  .1 , .12  , .1  1,  I4.MMM) 

OUTPUI  prints  out  both  the  Identification  headers  for  I. -table 
and  ELAMB- table  and  the  tabular  entries  themselves. 

OUTPUT  input  parameter  definitions,  origins,  and  units  are  given 
in  Table  29.  OUTPUT  output  parameters  are,  depending  on  the  value  of 
IK,  the  appropriate  heading,  and  either  the  state  or  co-state  entries 
for  each  year. 


Table  29 
OUTPUT  INPUT 


Pa rameter 

Def ini t ion 

Un  its 

Origin 

.IK 

St  a t <•/ co-st  at  e 

flag 

Dimensionless 

(integer) 

Ma  i n 

program 

JJ 

Number  of  entries 
in  the  state  and 

co-state 

Dimens  ionless 
( integer) 

Ma  i n 

program 

.12 

First  co-state 
entry  index 

Dimensionless 
( integer) 

Ma  i n 

program 

J3 

bast  co-state 
entry  index 

Dimens  ion  1 ess 
( integer) 

Main 

program 

MMM 

Output  flag 

Dimension  1 ess 
( integer ) 

Main 

program 

Subroutine  COSTS (T,M) 

After  optimal  bonuses  have  been  computed  and  the  C-table  con- 
structed as  a compact  representation  of  yearly  modified  states  (mod- 
ified by  adding  appropriate  linear  multiples  of  the  bonuses),  co-state 
and  bonuses,  COSTS  computes  the  yearly  bonus  cost,  the  penalty  cost, 
and  Hamiltonian  corresponding  to  the  optimal  solution.  It  also  com- 
pares total  system  costs  with  the  costs  accruing  to  the  evolution  of 
the  uncontrolled  systems  and  performs  the  net  gain  analysis. 
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COSTS  input  parameter  def  initions,  units,  and  origins  are  sliown 
in  Tab lo  JO.  COSTS  output  paramo tor  dot  init  ions  and  units  aro  given 
i n T ab  1 c ! 1 . 


Tab  1 o JO 
COST'S  INPUT 


i’a  rame  t e r 

Def init ion 

Un  i l s 

Origin 

T 

Problem  extent 

Years  (real) 

Ma  i n 

program 

M 

Problem  extent 

Years  (integer) 

Main 

program 

N 

St  a t e/ co-st a t e 

I)  i mens  i on  less 

Ma  i n 

program 

d intension 

(integer) 

N 1 

State  dimension 

Dimension  less 

Ma  i n 

program 

( integer) 

c; 

Mod i i i ed  slat e/ 

I Men , men/do  liar]' 

(IF  IX 

co-state  table 

( rea  1 ) 

uu 

Bonus  array 

Dollars  (real) 

Ma  i n 

program 

HI’ 

Supply  parameter 

Men/do  liar  ( rea 1 ) 

PROCO 

vector 


In  starting,  COSTS  prints  out  the  unmodified  C-table.  It  then 
calls  CFIX,  which  modifies  each  yearly  entry  of  the  C-table,  each  such 
entry  consisting  of  the  state/co-state  vector.  On  return  from  CK1X, 
COSTS  prints  out  the  modified  C-table.  COSTS  then  calls  STAR  for  the 
computation  of  costs  in  the  uncontrolled  system.  Next,  looping  over 
each  year  from  year  zero  to  year  M,  COSTS  computes  the  bonus  cost  as 
the  sum  of  squares  of  the  products  of  bonus  cost  multipliers  and 
bonuses.  The  penalty  cost  is  computed  as  a positive  semi -def in i to 
quadratic  form  in  the  difference  between  current  and  terminal  states. 

The  year's  total  cost  is  computed  as  the  sum  of  penalty  and  bonus  costs. 
The  gain  is  then  computed  as  the  difference  between  total  cost  in  the 
uncontrolled  and  optimally  controlled  system.  The  Hamiltonian  is  com- 
puted as  the  difference  between  the  vector  product  of  the  subsequent 
years  state  and  co-state  and  the  total  cost.  The  Hamiltonian  is  not 
computed  in  the  last  year  for  which  the  subsequent  year's  state 
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Table 11 

COSTS  OHTPU'I 

Par. irnet  er 

Del  in i t ion 

i'n  i t 

s 

TIMF. 

Current  year 

Year  (re 

al  ) 

liens  r 

Year's  bonus  cost 

I)ol  1 ars 

( rea  1 ) 

PCOST 

Year's  penalty  cost 

Do  1 1 ars 

( rea  1 ) 

TCOST 

Year's  total  cost 

Do  1 1 ars 

( rea  1 ) 

CAIN 

Year's  net  gain 
between  optimal 
system  and  uncon- 
trolled system  cost 

Do  1 1 a rs 

(real) 

11 

Year's  Hamiltonian 

Do  1 1 ars 

( rea  1 ) 

TOTALS 

Cumulative  total  for 
penalty  cost,  bonus 
cost,  total  cost,  un- 
controlled system 
cost  and  gain 

and  co-state  are  undefined.  Finally,  totals  for  penalty  cost,  bonus 
cost,  total  cost,  and  net  gain  are  accumulated. 

Subroutine  STAR(M) 

STAR  computes  the  yearly  states  and  associated  penalty  cost  for 
the  uncontrolled  system.  STAR  input  parameter  definitions,  units,  and 
origins  are  shown  in  Table  T2.  STAR's  output  is  the  G-STAR-table  of 
yearly  uncontrolled  manpower  levels  in  each  year  group  in  the  array 
of  associated  finally  cost  PCSTAR. 

STAR's  Initial  action  is  to  pre-set  uncontrolled  problem  dimen- 
sion. It  next  sets  up  the  state  transition  matrix  and  constant  input 
vector  and  initial  state  that  characterizes  the  uncontrolled  system. 
Then  it  sequentially  computes  each  year's  state  from  the  preceding. 
Having  constructed  the  G-STAR-table  of  uncontrolled  states  it  then 
computes  the  cost  quadratic  form  in  the  difference  between  each  year's 
state  and  the  desired  terminal  state. 
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Table 

STAR  1 

M 

NPUT 

Paramet  er 

Dt*f  in  i t inn 

Un  i t s 

Or  i g i n 

M 

Problem  extent 

Years  (integer) 

PROCO 

N 1 

Problem  state 

Dimens i on  less 

PROGO 

d linens  ion 

(integer) 

ALPHA 

St.ite  t ransi  t ion 

Dimens i on  1 ess 

PROCO 

paramet ers 

(real) 

Kl./.KRO 

Initial  state 

Men  (list) 

PROCO 

LI.T 

Final  state 

Men  (red) 

PROCO 
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Ap[H‘n<lix  A 


MATHKMAT I CAI,  TO  l>  ROC  RAM  SYMBOL  MAI’ 


Mat  In-mat  h a 1 Symbo  1 

DBMMP  Symbol 

1, 

t 

KI.VEC 

B 

t 

UU 

A 

A 

6 

BB 

C 

CC 

K 

S 

RR 

S 

F 

F2 

F 

F2 

Fl 

S' 

Fl 

ST 

S'K2  + BjFj 

PI 

(R  - S'F2S  - F*) 

W 

(R  - S ' F ,S  - F*)"1 

W 

p; 

PT 

A. 

A 1 

6’ 

BT 

B1 

BA 

Bi1 

BA 

p 

P2 

Q 

Q2 

Y 

Y 

H0BED1M  PAOI  BLAMUMOT  HIMD 
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M.ithem.it  ioa  1 Symbol 

DBMMP  Symbol 

7. 

BA 

I) 

B 

4> 

(AR.AI)  (complex) 

V 

H,  K1NV 

11  i 

II  I 

V-’ 

K 

<JiaR  ( Vi ! ) 

U1>0W 

0 - uT) 
diaK  - 

CVEC 

V diag  ViV1 

CTEMP 

1) 

C 

E11 

El  1 

F 

12 

El  2 

E21 

E22 

Lo 

EL ZERO 

lt 

ELT 

W1 

W 

W2 

TERM1 

W3 

BTZ 

(’t  , t 

G 

CB.t 

BC0ST 

CP.t 

PC0ST 

(T()T 

TC0ST 
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Appendix  B 

ECONOMICS  MODEL  LOCATIONS 


The  following  set  of  equations  is  used  bv  subroutine  COMPL'T  of 

the  DBMMT  to  compute  the  partials  matrix.  O'. . > . The  input  parameters 

I.  , i = 1 , 2,  m,  . u are  defined  and  described  in  Munch 

T , i 12 

( 1977). 


X1  " 


m / 2 

I L 


T,  i 


( B . 1 ) 


X 


2 


m 


l 

i=m/ 2+1 


L 


T,i 


„ - [v;» . 6,x;°  ] <■ 

Z = Qy 

p’  = p + 1 

q,  = 61(Q/X1)P' 
q 2 = 62(q/x2)p' 

W1  = w2ql/q2 

U*  - M - 1 


r * pZ/Q 
s = X2W2/(X1WJ) 


( B . 2 ) 

( B.  3) 

( B . 4 ) 

( B . 5 ) 

(B.6) 

(B.7) 

(B.8) 

(B.9) 

(B. 10) 

(B.  11) 


-ri7- 


I = Zt 
Z1  = r,,l 
'■2  = r'l  ) 

P = w2/z2 

Zj  j = z^ljCu’  - P's)/Q 
'■22  = z2cl)dJ'  " P'/s)/Q 

z,7  = z,q2(P  + P)/Q 

'll  = p [zu  ‘ zM/l] 

j22  = P [z22  “ z2/l  ] 

J12  = Plzi2  “ ziz2/l  1 

S.1 1 I . i = 1.  2,  m/2,  j = I,  2,  m/2 

•I , 7.  i = 1 . 2 m/2,  j = "'  + 1 , , m , 

‘ i ' \ 

j i = m/2+1 in,  j = 1,  2,  ...,  m/2 

Vl)7,  i = m/2+1,  ...,  m,  j = m/2+1,  ...,  m 


f H . 12) 
(B. I J) 

(B. 14) 

( B . 1 5 ) 

(B. IB) 

(B. 17) 

(B. 18) 

( B . 19) 

(B.20) 

( B . 2 1 ) 

( B . 2 2 ) 


The  economic  modi’ I and  model  parameters  underlying  these  equa- 
tions are  given  in  Munch  (1977). 


